clear
set more off

local yearlist $yearlist
local yearlistshort $yearlistshort
local varlist1 $varlist1
local varlist2 $varlist2

local svcyearlist $svcyearlist
local svcyearlistshort $svcyearlistshort
local svcvarlist1 $svcvarlist1
local svcvarlist2 $svcvarlist2 

*************************************************
*Individual Level Regressions
*************************************************

*Get data
use ${tempdir}/panelreg_inddata_30yr.dta, clear

*Variables for the relevant marriagemarket and childhood years
gen maryear = year2 - 10
gen chyear = year2 - 30

*Here, assign the state the match the panel controls as statefip, state of residence
*Change to bpl to match on birth place instead
gen fips = statefip

*Merge on the marriage market variables
gen year = maryear
sort fips year 
merge m:1 fips year using ${tempdir}/panel_controls.dta, keepusing(pcturban pcturb25 pctmale pctmale_wtyouth pctmale_wt1530 pctmfglab_ipums pctsvclab_ipums pctmfglabm_ipums pctmfglabf_ipums pctfarm_hhhd /*
	*/ pctfor pctforwt pctforwtf pctforwtm ginifarmsize_diff ginifarmsize_con pctrr15 pctrr40 logoccscore_nvrmar_2035 popperdwell valperfarm) keep(3) nogen

*Merge on the childhood variables - just schooling for now
replace year = chyear
sort fips year
merge m:1 fips year using ${tempdir}/panel_controls.dta, keepusing(pctsch_wt520 pctsch_wtm520 pctsch_wtf520) keep(3) nogen

*Merge on current outside logoccscore and rent
replace year = year2
sort fips year
merge m:1 fips year using ${tempdir}/panel_controls.dta, keepusing(logoccscore_nvrmar_3045) keep(3) nogen

replace pctrr15 = 0 if pctrr15 == .
replace pctrr40 = 0 if pctrr40 == .

*Create labels

label var pcturban "Urban"
label var pctmale "Percent Male"
label var pctmale_wtyouth "Pct Male Among White Youth"
label var pctsvclab_ipums "Pct Lab Force in Services"
label var pctmfglab_ipums "Pct Lab Force in Manufacturing"
label var pctmfglabm_ipums "Pct Male Lab Force in Mfg"
label var pctmfglabf_ipums "Pct Female Lab Force in Mfg"
label var pctfarm_hhhd "Pct Farmers Among HH Heads"
label var pctsch_wt520 "Pct of Whites 5-20 in School"
label var pctsch_wtm520 "Pct of White Males 5-20 in School"
label var pctsch_wtf520 "Pct of White Females 5-20 in School"
label var pctfor "Pct Foreign"
label var pctforwt "Pct Foreign Among Whites"
label var pctforwtf "Pct Foreign Among White Females"
label var pctforwtm "Pct Foreign Among White Males"
label var ginifarmsize_con "Gini of Farm Size"
label var pctrr15 "Land w/in 15 mi of RR"
label var pctrr40 "Land w/in 40 mi of RR"
label var popperdwell "Average Persons Per Dwelling"
label var valperfarm "Land Value per Farm"
label var logoccscore_nvrmar_2035 "Log Occscore Single Women 20-35"
label var logoccscore_nvrmar_3045 "Log Occscore Single Women 30-45"

foreach var in pcturban pcturb25 pctmale pctmale_wtyouth pctmale_wt1530 pctmfglab_ipums pctsvclab_ipums pctmfglabm_ipums pctmfglabf_ipums pctfarm_hhhd /*
	*/ pctfor pctforwt pctforwtf pctforwtm ginifarmsize_diff ginifarmsize_con pctrr15 pctrr40 pctsch_wt520 pctsch_wtm520 pctsch_wtf520 /*
	*/ logoccscore_nvrmar_2035 logoccscore_nvrmar_3045 popperdwell valperfarm {
	rename `var' `var'_stres
}

replace fips = bpl

*Merge on the marriage market variables
replace year = maryear
sort fips year 
merge m:1 fips year using ${tempdir}/panel_controls.dta, keepusing(pcturban pcturb25 pctmale pctmale_wtyouth pctmale_wt1530 pctmfglab_ipums pctsvclab_ipums pctmfglabm_ipums pctmfglabf_ipums pctfarm_hhhd /*
	*/ pctfor pctforwt pctforwtf pctforwtm ginifarmsize_diff ginifarmsize_con pctrr15 pctrr40 logoccscore_nvrmar_2035 popperdwell valperfarm ) keep(1 3) nogen


*Merge on the childhood variables - just schooling for now
replace year = chyear
sort fips year
merge m:1 fips year using ${tempdir}/panel_controls.dta, keepusing(pctsch_wt520 pctsch_wtm520 pctsch_wtf520) keep(1 3) nogen

*Merge on current outside logoccscore
replace year = year2
sort fips year
merge m:1 fips year using ${tempdir}/panel_controls.dta, keepusing(logoccscore_nvrmar_3045) keep(1 3) nogen

replace pctrr15 = 0 if pctrr15 == . & bpl < 100
replace pctrr40 = 0 if pctrr40 == . & bpl < 100

*Create labels

label var pcturban "Urban"
label var pctmale "Percent Male"
label var pctmale_wtyouth "Pct Male Among White Youth"
label var pctsvclab_ipums "Pct Lab Force in Services"
label var pctmfglab_ipums "Pct Lab Force in Manufacturing"
label var pctmfglabm_ipums "Pct Male Lab Force in Mfg"
label var pctmfglabf_ipums "Pct Female Lab Force in Mfg"
label var pctfarm_hhhd "Pct Farmers Among HH Heads"
label var pctsch_wt520 "Pct of Whites 5-20 in School"
label var pctsch_wtm520 "Pct of White Males 5-20 in School"
label var pctsch_wtf520 "Pct of White Females 5-20 in School"
label var pctfor "Pct Foreign"
label var pctforwt "Pct Foreign Among Whites"
label var pctforwtf "Pct Foreign Among White Females"
label var pctforwtm "Pct Foreign Among White Males"
label var ginifarmsize_con "Gini of Farm Size"
label var pctrr15 "Land w/in 15 mi of RR"
label var pctrr40 "Land w/in 40 mi of RR"
label var popperdwell "Average Persons Per Dwelling"
label var valperfarm "Land Value per Farm"
label var logoccscore_nvrmar_2035 "Lagged Log Occscore Single Women 20-35"
label var logoccscore_nvrmar_3045 "Current Log Occscore Single Women 30-45"

foreach var in pcturban pcturb25 pctmale pctmale_wtyouth pctmale_wt1530 pctmfglab_ipums pctsvclab_ipums pctmfglabm_ipums pctmfglabf_ipums pctfarm_hhhd /*
	*/ pctfor pctforwt pctforwtf pctforwtm ginifarmsize_diff ginifarmsize_con pctrr15 pctrr40 pctsch_wt520 pctsch_wtm520 pctsch_wtf520  /*
	*/ logoccscore_nvrmar_2035 logoccscore_nvrmar_3045 popperdwell valperfarm{
	rename `var' `var'_stbpl
}

* Regressions of ever married status
egen firstcode = group(sex first year)

gen int bpl2 = bpl
replace bpl2 = 9999 if bpl>100

gen byte foreign_born = bpl2>100 if bpl2~=.

sort year2 foreign_born

foreach var in pcturban pcturb25 pctmale pctmale_wtyouth pctmale_wt1530 pctmfglab_ipums pctsvclab_ipums pctmfglabm_ipums pctmfglabf_ipums pctfarm_hhhd /*
	*/ pctfor pctforwt pctforwtf pctforwtm ginifarmsize_diff ginifarmsize_con pctrr15 pctrr40 pctsch_wt520 pctsch_wtm520 pctsch_wtf520 /*
	*/ logoccscore_nvrmar_2035 logoccscore_nvrmar_3045 popperdwell valperfarm {
	by year2: egen hold = mean(`var'_stbpl)
	replace `var'_stbpl = hold if foreign_born == 1
	drop hold
	}

gen urb = urban==2 if urban!=.
label var urb "Urban dummy"


forvalues i = 1/4 {
	gen q`i' = qtl_score_wgt_full==`i'
	replace q`i'=. if qtl_score_wgt_full==.
}

foreach y in `yearlist'{
		gen y`y' = year2==`y'
}

forvalues i=1/4 {
	foreach y in `yearlist'{
		gen q`i'_`y' = q`i'*y`y'
	}
}

gen bpl1880 = bpl2 if year2 == 1880
replace bpl1880 = 7777 if year2 != 1880

local dep_ever_married = 1
local dep_logoccscore_dadinlaw_temp = 2

local sttype_stres = 1
local sttype_stbpl = 2

gen logoccscore_dadinlaw_temp = logoccscore_father_spouse 


cap postclose results
postfile results sttype depvar spec control qtile year b se using "$outputdir/graph_data.dta", replace

local c_stres ""
local c_stbpl "i.foreign_born##i.year2"

foreach state in /*stres*/ stbpl{
foreach var of varlist ever_married logoccscore_dadinlaw_temp {

* Ever married regressions -- no bpl controls
reg `var' `varlist1' i.age foreign_born if sex==2, cluster(firstcode)

forvalues i=2/4 {
	*local b`i'_1880 = _b[q`i']
	foreach y in `yearlist' {
		local b`i'_`y' = _b[q`i'_`y']
		local se`i'_`y' = _se[q`i'_`y']
	}
}

forvalues i=2/4 {
	foreach y in `yearlist' {
		post results (`sttype_`state'') (`dep_`var'') (1) (1) (`i') (`y') (`b`i'_`y'') (`se`i'_`y'')
	}
}

* Ever married regressions -- no bpl controls
reg `var' `varlist2' i.age foreign_born if sex==2, cluster(firstcode)

forvalues i=2/4 {
	local b`i'_1880 = _b[q`i']
		local se`i'_1880 = _se[q`i']
	foreach y in `yearlistshort' {
		local b`i'_`y' = _b[q`i'_`y']
		local se`i'_`y' = _se[q`i'_`y']
	}
}

forvalues i=2/4 {
	foreach y in `yearlist' {
		post results (`sttype_`state'') (`dep_`var'') (1) (2) (`i') (`y') (`b`i'_`y'') (`se`i'_`y'')
	}
}

* Ever married regressions -- full bpl controls
reg `var' `varlist1' i.age i.bpl2##i.year2 if sex==2, cluster(firstcode)

forvalues i=2/4 {
	*local b`i'_1880 = _b[q`i']
	foreach y in `yearlist' {
		local b`i'_`y' = _b[q`i'_`y']
		local se`i'_`y' = _se[q`i'_`y']
	}
}

forvalues i=2/4 {
	foreach y in `yearlist' {
		post results (`sttype_`state'') (`dep_`var'') (2) (1) (`i') (`y') (`b`i'_`y'') (`se`i'_`y'')
	}
}

* Ever married regressions -- full bpl controls
reg `var' `varlist2' i.age i.bpl2##i.year2 if sex==2, cluster(firstcode)

forvalues i=2/4 {
	local b`i'_1880 = _b[q`i']
	local se`i'_1880 = _se[q`i']
	foreach y in `yearlistshort' {
		local b`i'_`y' = _b[q`i'_`y']
		local se`i'_`y' = _se[q`i'_`y']
	}
}

forvalues i=2/4 {
	foreach y in `yearlist' {
		post results (`sttype_`state'') (`dep_`var'') (2) (2) (`i') (`y') (`b`i'_`y'') (`se`i'_`y'')
	}
}


local v3 pcturban
local v4 pctrr15
local v5 pctforwt
local v6 pctmale_wt1530
local v7 pctmfglab_ipums
local v8 pctsch_wt520
local v9 pctsvclab_ipums
local v10 logoccscore_nvrmar_2035
local v11 logoccscore_nvrmar_3045
local v12 popperdwell
local v13 valperfarm

foreach j of numlist 3/8 10/13{


reg `var' `varlist1' c.`v`j''_`state'#i.year2 i.age foreign_born `c_`state'' if sex==2, cluster(firstcode)

forvalues i=2/4 {
	*local b`i'_1880 = _b[q`i']
	foreach y in `yearlist' {
		local b`i'_`y' = _b[q`i'_`y']
		local se`i'_`y' = _se[q`i'_`y']
	}
}


foreach y in `yearlist' {
		local b99_`y' = _b[`y'.year2#c.`v`j''_`state']
		local se99_`y' = _se[`y'.year2#c.`v`j''_`state']
	}

forvalues i=2/4 {
	foreach y in `yearlist'  {
		post results (`sttype_`state'') (`dep_`var'') (`j') (1) (`i') (`y') (`b`i'_`y'') (`se`i'_`y'')
	}
}

foreach y in `yearlist' {
		post results (`sttype_`state'') (`dep_`var'') (`j') (1) (99) (`y') (`b99_`y'') (`se99_`y'')
	}


reg `var' `varlist2' c.`v`j''_`state'#i.year2 i.age foreign_born `c_`state'' if sex==2, cluster(firstcode)

forvalues i=2/4 {
	local b`i'_1880 = _b[q`i']
	local se`i'_1880 = _se[q`i']
	foreach y in `yearlistshort' {
		local b`i'_`y' = _b[q`i'_`y']
		local se`i'_`y' = _se[q`i'_`y']
	}
}

foreach y in `yearlist' {
		local b99_`y' = _b[`y'.year2#c.`v`j''_`state']
		local se99_`y' = _se[`y'.year2#c.`v`j''_`state']
	}

forvalues i=2/4 {
	foreach y in `yearlist' {
		post results (`sttype_`state'') (`dep_`var'') (`j') (2) (`i') (`y') (`b`i'_`y'') (`se`i'_`y'')
	}
}

foreach y in `yearlist' {
		post results (`sttype_`state'') (`dep_`var'') (`j') (2) (99) (`y') (`b99_`y'') (`se99_`y'')
	}

}	

forvalues j=9/9 {

/*
reg `var' `svcvarlist1' c.`v`j''_`state'#i.year2 i.age foreign_born `c_`state'' if sex==2, cluster(firstcode)

forvalues i=2/4 {
	*local b`i'_1880 = _b[q`i']
	foreach y in `svcyearlist' {
		local b`i'_`y' = _b[q`i'_`y']
		local se`i'_`y' = _se[q`i'_`y']
	}
}

foreach y in `svcyearlist' {
		local b99_`y' = _b[`y'.year2#c.`v`j''_`state']
		local se99_`y' = _se[`y'.year2#c.`v`j''_`state']
	}
	
forvalues i=2/4 {
	foreach y in `svcyearlist' {
		post results (`sttype_`state'') (`dep_`var'') (`j') (1) (`i') (`y') (`b`i'_`y'') (`se`i'_`y'')
	}
}

foreach y in `svcyearlist' {
		post results (`sttype_`state'') (`dep_`var'') (`j') (1) (99) (`y') (`b99_`y'') (`se99_`y'')
	}

	*/
	
reg `var' `svcvarlist2' c.`v`j''_`state'#i.year2 i.age foreign_born `c_`state'' if sex==2, cluster(firstcode)

forvalues i=2/4 {
	local b`i'_1880 = _b[q`i']
	local se`i'_1880 = _se[q`i']
	foreach y in `svcyearlistshort' {
		local b`i'_`y' = _b[q`i'_`y']
		local se`i'_`y' = _se[q`i'_`y']
	}
}

foreach y in `svcyearlist' {
		local b99_`y' = _b[`y'.year2#c.`v`j''_`state']
		local se99_`y' = _se[`y'.year2#c.`v`j''_`state']
	}
	
forvalues i=2/4 {
	foreach y in `svcyearlist' {
		post results (`sttype_`state'') (`dep_`var'') (`j') (2) (`i') (`y') (`b`i'_`y'') (`se`i'_`y'')
	}
}

foreach y in `svcyearlist' {
		post results (`sttype_`state'') (`dep_`var'') (`j') (2) (99) (`y') (`b99_`y'') (`se99_`y'')
	}
	

}	

}
}

postclose results

**** Make pictures *****

use "$outputdir/graph_data.dta", clear

label define state_lbl 1 "stres" 2 "stbpl"
label values sttype state_lbl

label define depvar_lbl 1 "ever married" 2 "dad in law's log occscore"
label values depvar devpar_lbl

label define spec_lbl 1 "no controls" 2 "bpl controls" 3 "pcturban" 4   "pctrr15"  5 "pctforwt"  6 "pctmale_wt1530" 7  "pctmfglab_ipums" 8 "pctsch_wt520"  9 "pctsvclab_ipums" 
label values spec spec_lbl

local state1 stres
local state2 stbpl

local l3 "% Urban Control"
local l4 "Railroad Access Control"
local l5 "% Foreign Control"
local l6 "% Male Control"
local l7 "% in Manufacturing Control"
local l8 "% Children in School Control"
local l9 "% in Services Control"
local l10 "Log Occscore Single Women Lagged Control"
local l11 "Log Occscore Single Women Current Control"
local l12 "Persons Per Dwelling"
local l13 "Land Value per Farm"

local l14 "Average Monthly Rent"
local l15 "Average House Value"


local stub1 evermarried
local stub2 dadinlawoccscore

local t1 "Probability of Ever Marrying"
local t2 "Spouse's Parental Log Occscore"

local y1 "-0.06(0.02)0.02"
local y2 "-0.01(0.01)0.03"

local yunit1 "Gap in Probability of Marriage"
local yunit2 "Gap in Parental Log Occscore"

gen birth_cohort = year-37.5

gen hi = b + 1.96*se
gen lo = b - 1.96*se


gen spec2 = .
replace spec2 = 1 if spec == 1
replace spec2 = 2 if spec == 8
replace spec2 = 3 if spec == 6
replace spec2 = 4 if spec == 13
replace spec2 = 5 if spec == 7
replace spec2 = 6 if spec == 12
replace spec2 = 7 if spec == 4
replace spec2 = 8 if spec == 3
replace spec2 = 9 if spec == 10
replace spec2 = 10 if spec == 5
replace spec2 = 11 if spec == 2

label define spec2labels 1 "Baseline"
label define spec2labels 2 "% in Sch", add
label define spec2labels 3 "% Male", add
label define spec2labels 4 "Val per Farm", add
label define spec2labels 5 "% Mnfg", add
label define spec2labels 6 "Ppl Per Dwell", add
label define spec2labels 7 "RR Access", add
label define spec2labels 8 "% Urban", add
label define spec2labels 9 "Sing Wmn Ear", add
label define spec2labels 10 "% Foreign", add
label define spec2labels 11 "Birth State", add
label values spec2 spec2labels 

	
foreach v in 1 2 {

cap drop birth_cohort_adj
gen birth_cohort_adj = birth_cohort - 1 if spec==1
replace birth_cohort_adj = birth_cohort if spec==2

graph twoway rcap hi lo birth_cohort_adj if qtile==4 & spec==1 & depvar==`v' & sttype==2 & control == 1, lcolor(gs4)  || /*
*/ rcap hi lo birth_cohort_adj if qtile==4 & spec==2 & depvar==`v' & sttype==2 & control == 1, lcolor(gs12)  || /*
*/ scatter b birth_cohort_adj if qtile==4 & spec==1 & depvar==`v' & sttype==2 & control == 1, msymbol(circle) mcolor(black) msize(medium) || /*
*/ scatter b birth_cohort_adj if qtile==4 & spec==2 & depvar==`v' & sttype==2 & control == 1, msymbol(diamond) mcolor(gs8) msize(medium) /*
*/ xtitle("Birth Cohort") ytitle("`yunit`v''") xlab(1840(20)1920) ylab(`y`v'') yline(0, lwidth(thin) lpattern(dash) lcolor(black)) /*
*/ title("Q4-Q1 Gap in `t`v''" /*"Women, 30-45"*/) /*subtitle("Birthplace Controls")*/ graphregion(fcolor(white)) /*
*/ legend(rows(1) order(3 4) label(3 "No Controls") label(4 "Birthplace Controls")) 

graph export "$outputdir/regression_`stub`v''_bplcontrol.pdf", replace
}

local v 1
local s 2
graph twoway rcap hi lo spec2 if qtile==4 & year == $endyear & depvar==`v' & sttype==`s' & control == 2 & spec == 1, lcolor(gs4) lwidth(medthick) || /*
*/ rcap hi lo spec2 if qtile==4 & year == $endyear & depvar==`v' & sttype==`s' & control == 2 & spec == 2, lcolor(gs4) lwidth(medthick)  || /*
*/ rcap hi lo spec2 if qtile==4 & year == $endyear & depvar==`v' & sttype==`s' & control == 2 & spec >= 3 & spec <= 13, lcolor(gs12) lwidth(medium) || /*
*/ scatter b spec2 if qtile==4 & year == $endyear & depvar==`v' & sttype==`s' & control == 2 & spec == 1, msymbol(circle) mcolor(black) msize(large) || /*
*/ scatter b spec2 if qtile==4 & year == $endyear & depvar==`v' & sttype==`s' & control == 2 & spec == 2, msymbol(diamond) mcolor(black) msize(large) || /*
*/ scatter b spec2 if qtile==4 & year == $endyear & depvar==`v' & sttype==`s' & control == 2 & spec >= 3 & spec <= 13, msymbol(square) mcolor(gs8) msize(medium) /*
*/ xtitle("Geographic Controls") xlab(1(1)11, val angle(45)) ytitle("`yunit`v''") ylab(`y`v'') yline(0, lwidth(thin) lpattern(dash) lcolor(black)) legend(off) /*
*/ title("Change in Q4-Q1 Gap in `t`v''") graphregion(fcolor(white))
graph export "$outputdir/Regression Controls Ever Marr.pdf", replace
	
local v 2
local s 2
graph twoway rcap hi lo spec2 if qtile==4 & year == $endyear & depvar==`v' & sttype==`s' & control == 2 & spec == 1, lcolor(gs4) lwidth(medthick) || /*
*/ rcap hi lo spec2 if qtile==4 & year == $endyear & depvar==`v' & sttype==`s' & control == 2 & spec == 2, lcolor(gs4) lwidth(medthick)  || /*
*/ rcap hi lo spec2 if qtile==4 & year == $endyear & depvar==`v' & sttype==`s' & control == 2 & spec >= 3 & spec <= 13, lcolor(gs12) lwidth(medium) || /*
*/ scatter b spec2 if qtile==4 & year == $endyear & depvar==`v' & sttype==`s' & control == 2 & spec == 1, msymbol(circle) mcolor(black) msize(large) || /*
*/ scatter b spec2 if qtile==4 & year == $endyear & depvar==`v' & sttype==`s' & control == 2 & spec == 2, msymbol(diamond) mcolor(black) msize(large) || /*
*/ scatter b spec2 if qtile==4 & year == $endyear & depvar==`v' & sttype==`s' & control == 2 & spec >= 3 & spec <= 13, msymbol(square) mcolor(gs8) msize(medium) /*
*/ xtitle("Geographic Controls") xlab(1(1)11, val angle(45)) ytitle("`yunit`v''") ylab(`y`v'') yline(0, lwidth(thin) lpattern(dash) lcolor(black)) legend(off) /*
*/ title("Change in Q4-Q1 Gap in `t`v''") graphregion(fcolor(white))
graph export "$outputdir/Regression Controls Husb Parental.pdf", replace
	
	
	